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SUMMARY 


A general treatment of parameter identification and stochastic control for 
use on ! alicopter dynamic systems is presented. The emphasis of the work 
reported herein is on rotor dynamic models, including specific applications to 
rotor blade flapping and the helicopter ground resonance problem. Dynamic 
systems which are governed by periodic coefficients as well as constant coeffi- 
cient models are addressed. The dynamic systems are modeled by linear state 
variable equations which are used in the identification and stochastic control 
formulation. The research presented addresses the pure identification problem 
as well as the stochastic ccntrol problem which includes combined identification 
and control for dynamic systems. The stochastic control problem includes the 
effect of parameter uncertainty on the solution and the concept of learning and 
how this is affected by the control's dual effect. The identification formula- 
tion requires algorithms suitable for on-line use and thus recursive identifi- 
cation algorithms are considered. The applications presented use the recursive 
extended Kalffl"”. filter for parameter identification which has excellent con- 
vergence for systems without process noise. Adaptive control results based 
upon the certainty equivalence principle are also presented for the rotor blade 
flap model. 


INTRODUCTION 

The concept of using active control for helicopter dynamic modes offers 
the potential for improved helicopter rotor operation. The alleviation of retreating 
blade stall, improvement in dynamic stability and reduction of vibration through 
active control are some of the potential improvements to be gained through 
active control. The current trend toward fly-by-wire computer systems and high 
frequency rotor control make the concept of active rotor control feasible. 

Higher harmonic swashplate control and individual blade control permit the use 
of active control techniques for stochastic control of dynamic modes. Active 
control of helicopter rotor systems requires the use of parameter identifica- 
tion of dynamic models and either fixed gain or adaptive control laws. The 
focus of this reasarch is on both pure identification of dynamic models and 
stochastic control as applied to helicopter dynamic systems. 



Until recently most of the research activity in parameter identification 
has focused on the problem of aerodynamic stability derivative extraction for 
handling qualities models [1,2, 3,4], The identification of stability and 
control derivatives for handling quality models is well understood, however 
by no means completely solved. Rotor blade ^denuif ication has also been done 
[5,6] with emphasis on wind tunnel applications. These rotor studies which 
address the pure identification problem have locused mainly on the physical 
model requirements and have used off-line type algorithms which do not formally 
account for random process noise excitation. Thus, applicability to general 
wind tunnel or flight test useage is limited. An on-line identification pro- 
cedure was presented in reference 7 for rotor blade mode identification from 
random wind tunnel turbulence excitation. A recursive maximum likelihood (RML) 
procedure was u. ed which is unbiased in the presence of both random measurement 
and process noise. Although the RML procedure was an on-line algorithm and 
suitable for random turbulence, ic is limited to single input-single output 
systems or autoregressive moving average (ARMA) models. 

The general state space identification problem is presented for helicopter 
rotor dynamic systems. Identification algorithms suitable for both random 
measurement and process noise are discussed for off-line and on-line use with 
emphasis on application to helicopter rotor dynamic models. 

The studies on control of helicopter systems have primarily focused on handling 

qualities autopilot design [8] with some research on the gust alleviation 
problem [9,10,11]. Recently, fully adaptive higher harmonic control (HHC) 
algorithms which include on-line identification and control have been developed 
and tested on simulations and in the wind tunnel [12,13,14,15]. These applica- 
tions are designed to minimize steady helicopter vibration which originate in 
the rotor. The acceptance and development by the rotorcraft community of high 
frequency rotor multicyclic hardware have made HHC suitable for flight vehicles. 

In addition, the present interest in helicopter individual blade control (IBC) 
[16,17,18] will further open up new areas of application for stochastic control 
and identification of rotor dynamic systems. Applications to rotor stall 
alleviation, improved rotor stability and gust alleviation are potential areas 
of application. 



This report presents the general state space formulation for helicopter 
rotor dynamic systems for parameter identification. The stochastic control 
formulation is presented which leads to fixed gain feedback controllers as well 
as adaptive control solutions. The stochastic control approach follows that 
presented in ref. 19 and 20 which includes the adaptive control properties: (1) 
caution, (2) probing, and (3) deterministic control. 

The theoretical background for identification and stochastic control is 
presented in Section II. The rotor blade flapping results which include 
periodic coefficients and reverse flow effects are presented in Section III. 
Parameter Identification results for the helicopter ground resonance model is 
presented in Section IV. 
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SYMBOLS 

Continuous time matrix used for ground resonance model defined 
in Eq. (3.24) 

time-varying matrix as a function of parameters 0 

scalar parameter and its mean value, respectively, used in 
Eq. (2.30) 

flapping model periodic coefficients defined in Eq. (3.6) - 
Eq. (3.8) 

value of flapping model periodic coefficients at y = 0 

time-varying control matrix as a function of parameters 0 

continuous time matrix used for ground resonance model defined 
in Eq. (3.25) 

measurement matrix (periodic) used for ground resonance model 
defined in Eq. (3.26) and Eq. (3.27) 

N-step cost function defined in Eq. (2.28) 

cost function evaluated at time step k 

denotes expected value 

denotes conditional expected value 

discrete transition matrix, approximated as in Eq. (2.6) 
general vector nonlinear function 

partial derivative of f w.r.t. z. (Jacobian matrix) 
nonlinear function associated with the state equation 
nonlinear function associated with the parameters, 0 
value of state transition matrix at time step k, (Eq. 2.58) 
longitudinal force excitation of support degrees of freedom, 

F 

§- , R = 20 ft 

m Q. R 
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lateral force excitation of support degrees of freedom, 
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discrete control matrix approximated as in Eq. (2.7) 
value of state transition matrix at time step k, (Eq. 2.58) 
measurement matrix 

measurement matrix associated with system state vector 

general nonlinear measurement function 

information set at time step k+1, defined in Eq. (2.38) 

time step number 

objective function (cost) 

optimal cost-to-go from time k to the end defined in Eq. (2.43) 
caution component of the tot.il cost 
deterministic component of the total cost 
probing component of the total cost 

backward solution of the matrix Riccatti equation computed by 
Eq. (2.59) 

time step number 

optimal CE control feedback gain computed by Eq. (2.57b) 

control gain for the deterministic one step example Eq. (2.34) 

stochastic control gain which includes caution given by Eq . (2.35) 
or Eq. (2.36) 

periodic coefficient associated with 8 of blade flapping given by 
Eq. (3.2), 1/sec 2 

periodic coefficient associated with g of blade flapping given by 
Eq. (3.3), 1/sec 

periodic coefficient associated with blade angle 0 of blade 
flapping given by Eq. (3.4) 

mass of the longitudinal support degree of freedom, (m x - 2767 slugs) 
mass of the lateral support degree of freedom (m^ = 2457 slugs) 
final time step 
number of state variables 

number of state variables plus parameters (n z * n x + ng) 
number of parameters 
covariance matrix 

conditional probability density function 
cost function weighting matrix on the state 
process noise covariance on the state equations 
process noise covariance on the parameters 
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cost function weighting matrix on the control (also used as noise 
covariance on the measurement) 

innovations covariance defined in Eq. (2,21) 

on-off switch to activate reverse flow terms in flapping equation 
continuous time index, sec 
integration time step, sec 
control vection 

deterministic and stochastic control solution defined in Eq. (2.50) 
through Eq. (2.52) 

continuous time random white process noise 
discrete time random white process noise sequence 

discrete time random white process noise sequence of x and 0 eqs 
mean value of v 

Kalman filter gain defined in Eq. (2.24) 
discrete time random white measurement noise 
state vector 

augmented state vector (x a * 0) 

updated state estimate at time k given measurements up to time k 
measurement vector 

state vector of states and parameters 

updated state estimate (states and parameters) at time k given 
measurements up to time k 

state error defined in Eq. (2.1b) 

flapping angle of blade, deg 
time step size 
Lock number, (y * 5) 
blade lag angle, deg 

cosing Fourier component of blade lag, deg 
sine Fourier component of blade lag, deg. 

longitudinal and lateral support damping ratio, respectively, 

Z critical 

blade inplane damping ratio, X critical 

rotor advance ratio, V/flR 

innovation sequence defined in Eq. (2.20) 

parameter vector 

rotor blade pitch angle, deg 
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correlation coefficient 
standard deviation 

standard deviation of a and b, respectively 
square root of the covariance between a and b 
covariance of x defined in Eq. (2.48) 

Jacobian matrix (Eq. (2.19)) 

general form of optimal control solution (Eq. (2.51)) 
azimuth angle of rotor blade, deg 

general form of control solution under the separation property 
(Eq. (2.52)) 

angular velocity of rotor, rad/sec (or RPM) 
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blade inplane natural frequency (t ft) 

longitudinal and lateral support natural frequency <* fl) 

denotes transpose 
denotes time derivative 
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THEORETICAL BACKGROUND 

This section presents the basic theoretical background for Identification 
and stochastic control of dynamic systems. A brief description of the identi- 
fication and stochastic control formulation is presented. This is followed by 
a more detailed treatment of the Identification problem. Finally, the funda- 
mental ideas of stochastic control theory are presented. 

Problem Description and Overall Approach 

The state variable representation is used to model helicopter rotor rotating 
coordinates and fixed axis coordinates. Rotating coordinates for either artic- 
ulated or hingeless rotors include blade flap, lag and torsion, which can also 
include rotor harmonic inflow coordinates. Fixed axis coordinates include 
vehicle body motion and landing support degrees of freedom for the ground 
resonance problem. A discussion of state variable models for articulated and 
hingeless rotors coupled to fixed axis degrees of freedom is given in ref. 1. 

The linearized state equations include periodic coefficients and can be repre- 
sented in either the rotating system or fixed axis system using the multiblade 
coordinate transformation. The linearized state equations will be used for the 
parameter identification problem and the stochastic control problem. 

The dynamic identification problem is based upon the linear state equations 
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i+1 


a 


i+i 


W *i + W u i + y i 
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( 2 . 1 ) 


where x i+i * s a 9tate vector 

Uj^ is a mxl control vector 

is a pxl vector of unknown parameters to be identified 
including periodic coefficients 

F i’^i are matr * ces as a function of the unknown parameter vector 0^ 
v^,w'^ are zero mean random white noise guassian sequences 
The measurement equation is 


vv 


x i + 
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where y^ is a rxl measurement vector 

is a rxn matrix of unknown parameters 0^ 

w^ is a zero mean white noise gauss; an sequence 

The identification objective is to determine estimates of the unknown para- 
meter vector 0£ and state x^ from the measurements y^; given control excitation 

V 

The stochastic control objective is to minimize with respect to the control 
a quadratic function of the state and control 


min E 


|(t *i *i + “I »i u i) + 4 q n j 


(2.3) 


where E{ } denotes expected value and 

Qi» are weighting matrices on the state and control, respectively 

The pure identification problem is concerned wi*-h Eq. (2.1) and (2.2). The 
pure control problem assumes the parameters are known In Eq. (2.1) and (2.2) 
and the optimal control u^ which minimizes Eq. (2.3) is to be designed. The 
stochastic control problem [19] simultaneously considers Eq. (2.1) and (2.3) 
and two such adaptive control designs are shown in block diagram form in fig. 1. 
Details of the parameter identification solution and stochastic control solu- 
tion are presented in the following sections. 


Parameter Identification Methods 


Parameter identification methods can be categorized as either an off- 
line or on-line method. Off-line methods usually require a number of passes 
or iterations over the data set being identified. On-line methods are often 
recursive and updated estimates are generated based upon past data and current 
measurements. Table 1 shows off-line and on-line identification methods which 
are asymptotically unbiased for systems contaminated by process noise, measure- 
ment noise, and both process and measurement noise. For systems with only process 
noise, the least square (LS) method yields identical results for either the off- 
line or recursive formulation. For systems with measurement noise only, the 
output error maximum likelihood method is suitable for off-line useage. The 
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recursive maximum likelihood (RML) method [20] and extended Kalman Filter (EKF) 
are both suitable for systems with measurement noise contamination. The RML 
method in [20] is based on single input-single output autoregressive moving 
average (ARMA) models while the EKF is formulated for state space models. For 
systems which include both process and measurement noise , the RML method or the 
"corrected" extended Kalman Filter [21] results in asymptotically unbiased para- 
meter estimates. 

Since the state space formulation is used throughout this analysis, the 
extended Kalman filter is used for parameter identification. Since the EKF 
results in biased parameter estimates in the presence of process noise, all 
simulation runs were contaminated with random measurement noise only. A brief 
description of the convergence characteristics of the EKF used for parameter 
identification is described below. 

The Extended Kalman Filter as an Identification Method and Its Convergence 
Properties . - The extended Kalman filter (EKF) is a well known method used for 
estimation of the state and parameters of a dynamic system. The popularity of 
the EKF as a parameter identification algorithm is due to the fact that the 
solution is recursive, is relatively simple to implement, and closely resembles 
the linear Kalman filter solution. It is known that the EKF often fails to 
converge to the correct parameter values under certain conditions. The first 
major proof of convergence of the EKF is as presented in Ref. 21 where modifi- 
cations to the EKF were presented which guarantees convergence to the true para- 
meter values. It was also established that the asymptotic convergence of this 
modified algorithm is as good as the off-line maximum likelihood (ML) algorithm. 
Since the off-line ML method is generally accepted as the "standard" for para- 
meter identification, (specially in aircraft applications), this result of Ref. 
21 is a significant justification for use of a "modified" extended Kalman 
filter. 

The algorithm being used for this research is the extended Kalman filter 
(without the modification of Ref. 21), since this algorithm was readily 
available and was previously implemented in the stochastic control algorithm 
used in this study. This algorithm provides the state and parameter estima- 
tion in the dual control algorithm intended for use in this research. The 
following summary is presented outlining the convergence properties of this 
algorithm. 

Case 1 - No process noise on the state model (e.g. no wind tunnel 
turbulence) 


/i 
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a) The EKF yields asymptotically biased estimates of the parameters 
if the parameter estimate is not close to the true value. 

b) The EKT yields asymptotically unbiased estimates of the parameters 
if the process noise covariance on the parameters is artificially 
set to a non-zero value. 

c) The covariance of the parameter estimate yields an accurate 
estimate of uncertainty of the parameter estimate for case lb. 

This estimate is a close approximation to the Cramer-Rao Error 
lower bound (ref. 34). 

Case 2 - with process noise on the state model 

a) The EKF yields asymptotically biased estimates, even when a non- 
zero value for the parameter process noise covariance is used. 

b) The modified EKF of ref. .’1 yields asymptotically unbiased 
estimates of the parameters and convergence is as good as off- 
line maximum likelihood r.dentif ication. 

The simulations performed in this research are done without process noise 
and thus the conclusions of case 1 aoove are applicable. The applications of 
the EKF algorithm to the flapping model and ground resonance model for the no 
process noise case is in agreement with the conclusions above for case 1. 

The identification application results have been excellent and the 
detail: ari presented in subsequent ructions. The EKF algorithm without the 
modification of ref. 21 as used in this effort is presented. A brief descrip- 
tion of the steps leading to the EKF and resulting equations are given in the 
next section. 

The Extended Kalman Filter Identification Method . - The use of the Extended 
Kalman Filter for the simultaneous state estimation and parameter identification 
for a linear system with unknown parameters based on noisy output observations 
is presented. While this is not an optimal algorithm, its recursive nature and 
ease of implementation, especially in view of the procedure presented here, make 
it appealing. Other algorithms, like the Maximum Likelihood (ref 22 and ref. 23) 
are known to be superior and since they are off-line algorithms are not 
Investigated here. 
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Consider the linear system described in continuous time by the equation 
for the n^ dimensional state vector 

x(t) - A(0(t))x(t) + B(0(t))u(t) + v ( t ) (2.4) 

where u is the n^ diemsnional control, v the process noise and 0 represents the 
ng-vector of unknown parameters. The discretized state equation with sampling 
interval A, sufficiently small, i.e., an order of magnitude below the smallest 
time constant or natural period (inverse of natural frequency) . is 

x(k+l) * F(0(k))x(k) + G(0(k))u(k) + v x (k) (2.5) 

where 0(k) is the value of the parameter vector during period k, 

F(0(k) ) = I + A(0(k))A (2.6) 

G(0(k) ) - B(0(k))A (2.7) 

and the control u is assumed to have the constant value u (k) during the k-th 
sampling interval. The discrete-time process noise v x (k) is assumed to be a 
zero-mean uncorrelated sequence with known variance matrix Q x . 

Let the measurement equation be 

y(k) - H X (k)x(k) + w (k) (2.8) 

where H x (k) is the, possibly time-varying, measurement matrix and w(k) the 
measurement noise. The measuiement noise is assumed to be a zero-mean uncorre- 
lated sequence with known covariance matrix R. 

The problem is to estimate simultaneously the state x and the parameter 
vector 0. If the parameters were known the problem would be a linear one and the 
Kalman Filter provides the linear minimum mean square error estimate. However, 
when 0 is not perfectly known a conmon approach is to use the "Extended Kalman 
Filter" which estimates the "augmented" state 

z(k) * 


x (k) 
0(k) 


(2.9) 
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n * n + n 0 (2.10) 
z x 0 

Since elements of the augmented state vector z, namely the components of 
the proper state x and components of the parameter vector multiply each other 
in Eq. v.2.5) it is obvious that one faces a nonlinear problem. The approach 
via linearization is described in the next section. 


The evolution of the augmented state can be written as 



x(k+l) 


F (0(k) ) 0 


x(k) 


G(0(k)) 


v X (k) 

z(k+l) = 

0(k+l)_ 

SS 

* X 


0(k! 

+ 

0 _ 

u(k) + 

v 0 OO_ 


where the parameter vector is modelled as a Wiener process driven by the zero- 

0 0 

mean uncorrelated noise sequence v (k) with an assumed covariance matrix Q . 

If the parameters are assumed constant then Q® ** 0. 

Eq. (2.11) can be written as 

z (k+1) = f[z(k), u(k) ] 4 v(k) (2.12) 

and the measurement equation (2.8) becomes, in terms ot the augmented state z 

y (k) = H(k) z(k) + w(k) (2.13) 

The Extended Kalman Filter equations are obtained as follows. Assume at 
time k one has the estimate z(k|k) and the associated covariance matrix P(k|k). 
Then the predicted stare is obtained by the following first order expansion 

z(k+l' * f [z(k|k) ,u(k) ] 4 - f z [z(k|k) ,u(k)] [z(k) - z (k j k) ] + v(k) (2.14) 

where f is the Jacobian of the function f with respect to z, evaluated at the 
latest estimate z(k|k). The evaluation of the Jacobian and an automated computer 
implementation of it are discussed in the Appendix. 

The predicted state z(k+l|k) is obtained by taking the expected value of 
(2.14) conditioned on the observations up to and including time k, yielding 

z(k+l|k) = f[z(k|k),u(k)] (2.15) 


where the error 
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z(k|k) • z(k) - z(k|k) (2.16) 

is assumed zero-meuti and thus does not appear In (2.15). 

Subtracting (2.15) from (2.14) yields 

z(k+l|k) A z(k+l) - z(k+l|k) ■ f z [z(k|k), u(k)|lz(k) - z(k|k)] +v(k)(2.l7) 

The covariance of the p^ec’c.ied state is thus 

P(k+l|k) - E[z(k+1 | k) z' (k+1 |k) ] - $(k)P(k|k)*'(k) +Q (2.18) 

where 

*(k) $ f [z(k|k), u (k)l (2.19) 

z 

The covariance propagation equation (2.18) is the same as for linear systems 
with the Jacobian taking the place of the transition matrix. 

The covariance of the Innovation 

v(k+l) & y (k+1) - H(k+l)z(k+l|k) (2.?0) 


is 


S(k+1) - H(k+l)P(k+l|k)H' (k+1) + R (2.21) 

as in the linear case. 

The state and covariance update are given by the ..-andard equations 


z(k+l j k+1) - z(k+l|k) + W(k+1)\'<.K+1) (2.22) 

P(k+1 j k+1) - [I - W(k+l)H(s-rl)]P(k+l|k) (2.23) 

where the filter gain is 

W(k+1) - P(k+l|k)H' (k+r=3~ i (k+l) * P(k+1 | k+l)H* (k+l)R _1 (2.24) 

The alternate expression for the updated covariance is numerically more stable 
P(k+1 1 k+1 ) - P(k+1 j k) - P(k+i|k)H’(k+l)s" i (k+l)H(k+l)P(k+l|k) (2.25) 

than (2.23). 


A similar linearization can be carried out for the measurement equation 
when it is nonlinear. 
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Stochastic Control Theory 

The problem of stochastic concrol - control of uncertain systems where the 
uncertainty is modeled in probabilistic terms - consists, at first sight, of 
the following subproblems 

(i) extraction of information from the system (estimation/identif ication) 

(ii) control of the system based on this information 

The topic of this section is the intimate connection between the above two 
functions, which, in general, cannot be separated. 

The next section introduces the assumptions of the Bayesian Stochastic 
Control framework.. The effect of parameter uncertainty on the solution of a 
simple stochastic control problem, which usually leads to "caution" on the part 
of the controller, is shown to be not necessarily intuitively predictable. This 
section also introduces the concept of learning and how this is affected by the 
control's dual effect. The Principle of Optimality for stochastic systems 
together with the resulting stochastic dynamic programming equation are presented. 
As a consequence of the Principle of Optimality for stochastic systems the 
optimal controller has to anticipate (causally) the future system uncertainties. 
Specifically, the availability of future observations and the accuracy of future 
state/parameter estimates enter into the controller's decision about the current 
control value. This opens the way for actively adaptive controllers that can 
reduce system uncertainties via the dual effect and then enhance state/parameter 
estimation to ultimately improve control performance. Finally, the intimate 
connection between the control's dual effect and Certainty Equivalence property 
of the optimal stochastic control for a class of systems is given. 

The Basic Modeling Assumptions in Stochastic Control - The Bayesian 
Approach . - In the Bayesian approach for the control of uncertain systems the 
uncertainties are modeled in probabilistic form. Imperfectly known initial 
states and system parameters are modeled as random variables while input and 
output disturbances are modeled as white noise sequences. The distinction 
between states and parameters is somewhat fuzzy - the latter are varying much 
slower than the former. In view of this they are sometimes lumped together 
in what is called augmented state. 

The general model of a discrete-time system with unknown (and possibly 
time-verying) parameters can be written as 


14 



OWGIN/U. PAGE IS 
0F POOR QUALITY 

x(k+l) = f x [k, x(k), 0(k), u(k), v X (k)] (2.26a) 

0(k+l) - f°[k, 0(k), v°(k)] (2.26b) 

y(k) - h [k, x(k) , w(k) ] (2.26c) 


wb re x is the (proper) state vector, 0 the parameter vector, u the control, v x 
c e state process noise, v G the parameter process noise, y the measurement (out- 
pi t) and w the measurement noise. With the augmented state x - r v * 
be written as 


[x 1 0']\ (2.26) 


x a (k+l) - f [k, x a (k), u(k), v(k) ] (2.27a) 

y(k) = g [k, x a (k), w(k) ] (2.27b) 

The reason for the noises being modeled as white is to have x a Markov process 

a i k k- 1 

and then (and only then) its pdf (probability density function) p[x (k) | Y ,U ] 

conditioned on the available measurements Y = {^(j ) , j=0, . . . ,k} and past con- 
trols U* 1 ■ {>i(j) ,j=°0, . . . ,k-l} summarizes the past and constitutes the "infor- 
mation state" for the stochastic problem (see Ref. 24). 

The Bayesian approach consists of minimizing the expected value of a cost 
function (or maximizing a utility function). While other approaches like minimax 
or worst distribution have been also considered, they are less tractable for 
closed form solution. Consider, for system (2.26), a fixed end-time problem of 
N steps with uhe usual cost 

N-l 

C(0, X N , U N_1 ) = c [N, x(N) ] + £ c[k, x(k) , u(k)] (2.28) 

k-0 

N 

where X » {x(k) , k=0 , . . . ,N-l} . The stochastic control problem is then 


min E [ C ] (2.29) 

„N-1 

T'i order for the above expectation to exist it is necessary that every 
variable entering directly or indirectly into (2.28) be either deterministic 
'i.e., perfectly known) or random with a suitable pdf attached to them. Thus 
"unknown constant parameters" are modeled as a realization of random variables 
according to their pdf. 

Tu thermore the implications of the expectation in (2.29) are that we want 
to find the optimal control policy over the ensemble of all possible initial 
3 tates parameters and disturbances. The fact that the optimum is sought for 
the ensemble of possible parameters (i.e. different realizations of the plant) 


is an important point. This has an implication on the control's "adaptation" as 
will be discussed later. 


Parameter Uncertainty and Learning .- Representation (2.26) of a system is 
more convenient than the augmented system (2.27) when (2.26a) is linear in x, 
u and v x , i.e., with known parameters the system would be linear. Consider such 
a system (scalar for the purpose of this discussion) 

x(k+l) * a x(k) + b u(k) + v(k) 


0(k+l) * [a b] * - Q(k) 
y(0) - x(0) 

where v(k) is a zero-mean white noise. 

First let N * 1, i.e., a single step problem with cost 

C = x 2 (l) 

The optimal stochastic control for this (static) problem is 

u S (0) * -L S x(0) 

where the control gain for the stochastic situation is 


(2.30) 


(2.31) 


(2.32) 


L S . E[a b 1 
E[b 2 ] 


a b * p ab a a °b 

-2 2 
b + O, 

D 


(2.33) 


If c ■ ct. ■ 0 the control gain for the corresponding deterministic situation is 
a d 


t D a 
L “ — 


(2.34) 


The effect of uncertainty can be a decrease or an increase in the control 

gain: 

(i) If o - 0 
a 


.S a b , .D 

b + CF 

D 


(2.35) 


i.e., "caution" will be exercised by the controller because of parameter 
uncertainty. 

(ii) If a ■ b ■ 1, a b”l 
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(2.36) 


S D 

and clearly one can have L > L 

Thus the effect of the parameter uncertainty has to be analyzed individually 
for each problem - even in a simple static case the results are not necessarily 
intuitively predictable. 

In multistage control problems for systems with unknown constant parameters 
the initial parameter uncertainty is modeled by a prior pdf p[0|l^] where 1^ is 
the initial information. The initial control u(0) will (among other things) 
account for the fact that it is applied to a system with parameter 0 "drawn'* 
from the prior pdf. Having assumed the parameter as time-invariant one can 
expect that the initial uncertainty about the parameter's true value will be 
reduced in the course of the control process - the controller can "learn" it. 
Thus, as new information is gathered from measurements, the controller can adapt 
itself to the particular system it is controlling. 

Consider again system (2.30) in the multistage situation with known a and 
unknown b. Assume the pdf of b to be Gaussian with mean b^ and variance 
and the variance of v to be Due to the linearity of the estimation problem 

for b, one obtains the recursion for its variance 


2 r i k-fl , 2 r 2,. . 2 -2 , -1 

° k+1 * var[b|I ] « o [u (k) + a a k ] 


(2.37) 


where 


r k+l 


[Y k+1 , U k } 


(2.38) 


is the information set: ac k+1. 

Thus the prior uncertainty is reduced but it also depends on the control - 
the control has a dual effect (Ref. 25); 

(i) has an effect on the state 

(ii) affects subsequent information accuracy - parameter identification in 
this case. 

This suggests that the control could be used to enhance the identification 
process while controlling the system by "probing" it (Ref. 25) in order 
to increase the accuracy of subsequent control actions. Such a controller is 
called dual or actively adaptive . 
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The Principle of Optimality for Stochastic Systems . - The Principle of 
Optimality (Ref. 26) for the control of a stochastic system can be stated 

as follows: at any time, whatever the present information and past decisions, 
the remaining decisions must constitute an optimal policy with regard to the 
current information set. 

Since the principle of optimality states that every end part of the decision 
process must be optimal, the multistage optimization has to be started from the 
last stage. The last decision, u(N-l) , must be optimal with regard to the in 4 ,r- 
mation set available when it has to he computed, i.e., it will be obtained from 
the functional minimization 

min E(c| I N_1 ) (2.39) 

u(N-l) 


where C is the cost for the entire problem. 

The next to the last decision, u(N-2) 

N-2 

1) must be optimal with respect to (w.r.t.) I and 

2) is to be made knowning that the remaining decision u(N-l) will be 
optimal w.r.t. I ^1 

Thus, the (functional) minimization that yields the decision function at 
N-2 is 


min E min E(C | I N_1 ) | £ N “ 2 
u(n-2) u(N-l) 


(2.40) 


and it uses the result of the functional minimization (2.39). 

Note that the outside averaging in (2.40) is over y(N-l) using the condi- 
tional density 


p[y(N-l) |l N_2 , u(N-2) ] 


(2.41) 


parameterized by the control at N-2. Since this measurement is not yet 
available when u(N-2) is to be computed but it will be available for u(N-l) 
it is "averaged out" in (2.40). 

The above-described last two steps are entirely similar to the "preposterior 
analysis" technique from the operations research literature discussed, e.g., in 
reference 27. This technique is usually formulated in the 

following context. The first decision (here u(N-2)] is for information gathering 
by an experiment from which a posterior information will result [here y(N-l)] 
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that will be used to make the la3t decision [here u(N-l}]. The prior (to the 
experiment) probability density of the (posterior) result of the experiment is 
called the "preposterior density" and in the present probl^an this is (2.41). 
Thus, one can say that preposterior analysis, which is "anticipation" (in a 
statistical sense, i.e., causal) of future information is a consequence of the 
principle of optimality. 

This "causal anticipation" of future information is the key issue in 
obtaining a dual controller. 

The extension of (2,41) to the full N-stage process yields the optimal 
expected cost starting from the initial time as 


* 0 

J (0,1 ) - min 
u(0) 



min E min E(C|l N_1 ) |l N_2 
u(N-2) u(N-l) 



(2.42) 


where 1^ is the initial information. Note that this equation does not assume 
any particular form for the cost function ?. 

For the additive cost given by (2 28' one obtains from (2.42), after some 
manipulations (see, e.g., Ref. 24), the backward recursion known as 

Bellman's equation or stochastic dynamic programming 


J (k,I k ) ■ min E (c[k,x(k) ,u(k) ] 
u(k) 


+ J*(k+1, I k+1 ) |l k } k*N-l, . . . ,0 (2.43) 

^ lc 

where J (k,I ) is the optimal cost-to-go from time k to the end and its 
dependence on the available information set at k is explicitly pointed out. 

The terminal condition for (2.43) is 

J*(N,I N ) - E {c[N,x(N)] |l N ) (2.44) 


where the last measurement is irrelevant since it is averaged out immediately. 

The deterministic dynamic programming equation - similar to (2.43) except 

k. 

without the expectation and with I replaced by x(k) - can be solved only when 
an explicit expression of the optimal cost-to-go can be obtained recursively. 
Otherwise numerical techniques have to be used, but they are of limited use- 
fulness due to the "curse of dimensionality" - the number of required quantiza- 
tion points in the state space increases exponentially with the dimension of 
the problem. 
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The stochastic dynamic programming equation (2.43) has the additional 
problem of averaging over the next measurement [as indicated in (2.41)] and 
over the state conditioned on the current information set . All this causes 
the practical usefulness of (2.43) to be of limited value and indicates the 
need to find suboptimal controllers by suitable approximation of (2.43). The 
discussion from the last section indicated the two features that a stochastic 
controller should have 

(i) caution 

(ii) probing via the dual effect 

The approximations of the stochastic dynamic programming fall in the 
following two classes. 

1) Feedback Type A ZgotUthito: In this case the control depends only on 

the current information 

u(k) “ u(k, I k ) (2.45) 

but does not use the prior statistical description of the future 
posterior information 

p[y(j+l)|l j ], j >k. (2.46) 

2) Closed- Loop Type AZgoxZthm : Such a controller utilizes feedback 
(2.45) and anticipates future feedback via (2.46), i.e., that the 
loop will stay closed. 

It is clear that an approximation of (2.43) that has both the features of 
caution and probing will have to be of the closed-loop type. 

The Control’s Dual Effect and Certainty Equivalence . - Based on the discussion 
of the previous sections, the control is said to have a dual effect if the 
uncertainty about the stochastic system's state (augmented to include the unknown 
system parameters; superscript a is deleted for simplicity) depends on past 
control values. This can be formalized as follows. Let 

x(k|k) - E[x(k)|l k ] (2.47) 

Z(k|k> - E{[x(k) - £(kjk)][x(k) - «(kjk)]’|l k } (2.48) 

k-1 

Then, if the covariance (2.48) of the augmented state does not depend on U 
the control has no dual effect (or second order) - the control is neutral. 
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Otherwise it has the dual effect since it affects the uncertainty about the 
augmented state x(k). 

It should be pointed out that the dual effect of the control is not limited 
to systems with unknown parameters but occurs in general in nonlinear systems. 

The dual effect of the control and the Certainty Equivalence (CE) property 
of the optimal stochastic control have been shown to be inter-related for a 
class of problems (Ref. 28). 

The CE property is defined as follows. Consider the stochastic problem 
(2.29J with an arbitrary cost for system (2.27). Let the corresponding deter- 
ministic problem be 


min C (2.48) 

UN-1 


for system 


x(k+l) - f[k, x(k) , u(k) , v(k) ] 
y(k) - x(k) 


(2.49) 


where the noise has been replaced by its mean v. The solution to this deter- 
ministic problem in feedback form can be written as 

u°(k) - <J>[k, x(k)J (2.50) 


If the optimal solution to the stochastic problem (2.4), (2.2) is 


u S (k) - <Hk, x(k|k) ] 


(2.51) 


l.e., it is the same as the deterministic one except that the state is replaced 
by its conditional mean, then the stochastic controller has the CE property. 

In other words the "control law" $ is the same, only the state estimate replaces 
the state. As will be indicated later the CE property holds under rather 
restrictive conditions. 

If the control law is different, but still uses as input only the state 
estimate 

u S (k) - Mk, x(k j k) ] (2.52) 

one has the weaker property of separation - there is an estimator cascaded with 
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a controller. A problem where separation holds but not CE is the so-called 
Linear-Exponential-Quadratic-Gaussian (LEQG) problem (Ref. 29). 

The following result (Ref. 28) connects the control's dual effect to the 
CE property. 

Theorem . For the linear system with additive zero-mean white noise 


x(k+l) - F(k) x(k) + G(k) u(k) + v(k) 
with observation of the general form 

y (k) - h[k, x(k) , w(k) ] 


and quadratic cost 


N-l 

C - x' (N)Q(N)x(N) x* (k)Q(k)x(k) + u' (k)R(k)u(k) 

k-0 


(2.53a) 


(2.53b) 


(2.54) 


the optimal stochastic controller that minimizes the expected value of (2.54) 
has the CE property if and only if the control has no dual effect (of second 
order) . Then 

u S (k) - L(k) x(k|k) (2.55) 


where L(k) is the same feedback gain as in the corresponding deterministic 
problem. 

Note that., in general, when there are nonlinear observations the estimation 
accuracy £(k|k) is control-dependent. In the case of linear observations, how- 
ever, there is no dual effect and (2.55) holds. This is the reason for the 
popularity of the LQG problem whose solution Is the cascaded linear deterministic 
feedback with the Kalman filter that yields x(k|k). 

The CE result has been used in the context of adaptive control as follows. 
For the system with unknown parameters 


x(k+l) - F(k,0) x(k) + G(k,0) u(k) + v(k) 
y(k) * H(k) x(k) + w(k) 


(2.56) 


with quadratic cost (2.54) one can break down the problem into two parts 

(i) estimation/ identification 

(ii) control (using latest estimates) 
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The controller can be taken as 

u(k) - L[k t 0(k)] x(k| k) (2.57a) 


where G(k) is the estimate of 0 at time k and it is used in the feedback gain L 
as if it were the true value. This control, can be called as Heuristic Certainty 
Equivalence (HCE) because the CE property does not hold in this problem. Note 
that this controller does not account for the uncertainty in § (no caution) nor 
does it account for the opportunity of probing the system. Nevertheless, con- 
trollers of this type, which include the Self-Tuning Regulator have been used 
very successfully (Ref. 30). 

The optimal feedback gain for the HCE controller is 

Sc • (R k + C k *fcfl V" 1 C k bc+1 F k (2 - 57b > 

where 


L k - L(k, 0(k)) 

F k - F(k, 0(k)) (2.58) 

G k - G(k, 0(k) ) 

The matrix is given by the backwards iteration solution of the Rii-catti 
equation 

■Sc ■ f ; ivi - ■'k+i k + c ; Vi v 1 s F k + <>k 

S'On (2 - 59) 

The optimal feedback gain given by Eq. (2.57) and (2.59) is the same as in the 
deterministic problem. 

The structure of the solution for the optimal control u(k) of Eq. (2.57) 
is the same as the LQG solution (Ref. 31). 

The HCE type controller is of the feedback or passively adaptive type. A 
closed-loop type "dual" controller, which is actively adaptive, was developed 
in 1'se and Bar-Shalom (Ref. 32, 33, and 19). For further discussion on the 
quantification of the caution and probing concepts see reference 24. 
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ROTOR BLADE FLAP IDENTIFICATION AND CONTROL 

l^s section summarizes the results of application of the EKF identification 
method and the stochastic control solutions as discussed in the previous section. 
The methods are applied to a single flapping rotor blade modeled by both 
constant and periodic coefficients. 

Rotor Blade Flap Model 

The rotor blade flap model used in the simulation analysis includes periodic 
coefficients and the effects of reverse flow for rotor azimuth ip greater than 180° 
and less than 360°. The equations are linear in the state and deled explicitly 
as a function of advance ratio in order to investigate the Incre. ed effects of 
periodicity and reverse flow as advance ratio is increased. 

e rotor blade flap model equations were taken from ref. 38 and are shown 
ir Eq . (3.1) through Eq. (3.4). 

6 + A - y[M^B + MgS + m^ 0 r ] (3.1) 

M 0 - -y cosiKl/6 + 1/4 y sin^ -s(ip) sin\) (2 2 (3.2) 

P o 

y 4 4 

Mg ■ (- 1/8 - 1/6 y sin^ - s(^) s * n 40 ^ (3.3) 

4 

Mg - (1/8 + 1/3 y sinlp + 1/4 y 2 sin 2 if/ - s(ip) ^ sin%) fi 2 (3.4) 

R 

The reverse flow terms are effective for 180° < ij,' < 360° using the on-off 
switch 


s(p) 



0 < <p <_ 180° 
180° < t p < 360° 


The blade flap model was simulated at ft » 320 RPM with a Lock number y * 5. 
Advance ratios of y ■ 0 to y ■ .8 are investigated. 

The discrete state space representation of Eq. (3.1) through Eq. (3.4) is 


x 1 (k+l) 

m 

1 

L w 


x 1 (kT 

+ 

0 

x 2 (k+l) 


*21 

a 22 


x 2 (k) 


*21 


u(k) 


(3.5) 
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where the discrete coefficients , *22 ani * *’21 are using the first 

order continuous to discrete approximation of Eq. (2.6) and (.2.7). 

The A and B matrix elements are 

3 

a 21 " *21 ~ c08 ^ “ ^3 cosd. + G^ 8 (||)) sin H 1 coatj' (3.t>) 

a 2 ? " *22 ^2 8inl ^ " O 5 8 (H)) sin% (3.7) 

b,i - sinH) + 0^ sin^H 1 - s(H0 Og sin*<|> (3.8) 

where the reverse flow switch s(l()) is defined as before. 

The parameter vector for identification is 

0 * (* 2 i *22 ^21 ^1 ^2 ^3 ^4 ^5 ^6 ^7 ^ 8 ^ (3.9) 

where p. rameters 0. through 0 U are coefficients of oeriodlc terms. 

1 o 

In order to assess the importance of the individual periodic terms the 
magnitude of the parameter vector 0 is tabulated in Table 2 for increasing 
advance ratio from 0 < y < . 8 . For p ■ 0 there is no periodicity, p • .2 shows 
weak periodicity and negligible reverse flow, p ■ .5 shows strong periodicity 
and some reverse flow, s~d p * .8 shows strong periodicity and reverse flow, 

Thi next section summarizes the results obtained for pure identification 
over an advance ratio range of p • 0 to p • . 8 . Stochastic control applications 
are then presented. 

Parameter Identification Results 

The parameter identification formulation was discussed in section II which 

presented the state variable description (Eq. (2.1) and Eq . (2.2)) tor use with 

the EKF. The state variable flapping model used in the identification is given 

by Eq. (3.3) with coefficients defied In Eq . (3.6) to Eq . (3,8). 

All pure identification results assume a pulse input excitation 0^ for a 

duration of one sample period (At * -005). The magnitude of the input is such 

as to vield an initial condition on • 33 for adequate system excitation. 

sec 

This corresponds to an input amplitude of 0 * 10 deg. lhe time history response 

* 

for and 8 is shown in figure 2. The single blade measurement (3 Is used in all 
runs. The damping ratio is C * .316 (y/lb) and for SI ■ 33.51 rad/sec (320 RPM) , 
the damping factor is 10.47. The flapping time constant is .0933 seconds. 
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Identification results at hover, y » 0 . - Since there is no periodicity 
for y »■ 0, only two parameters are identified from the initial condition response 
of figure 2. The parameter vector is 

0 T » [a 21 a 22 l (3.10) 

Two noise levels were used in the simulation, a high noise case 2o y «- .08 deg 
and a low noise case 2o y ■ .008 deg. The high noise case is approximately a 
20% noise-to-signal ratio and the low noise case is approximately 2% noise on the 
flapping signal 6 . 

Identified parameter convergence is shown in figure 3 for the high noise 
case (20% noise) and in figure 4 for the low noise case (2% noise). In both 
cases convergence to the simulation value occurs in less than 1/2 of a rotor 
revolution. Parameter convergence occurs in a 1/4 of a revolution for the low 
noise case. Initial starting values for the parameters a^ and a 22 and standard 
deviations are also shown. The standard deviations are computed from the Cramer- 
Sao lower bound, which for no process noise is accurately represented by the 
covariance equation of the Kalman Filter (ref. 34) . 

Identification results at |i 3 .2 . - The flapping model of Eq. (3.6) through 
Eq. (3.8) shows analytic expressions for the state variable A-matrix elements 
and a 2 2* These elements are plotted in fig. 5 for advance ratio y ■ .2. Small 
periodicity is shown vs rotor revolution which also show negligible reverse flow 
effects. Since, the reverse flow effects are negligible for this case, only the 
parameter associated with sin^ cosijj is identified. The parameter vector for 
identification is 


0 T * [Qj. 0 2 0 3 ] (3.11) 

The constant terns a 2 ^ and a 22 in Eq. (3.9) are assumed known, being identified 
previously at y = 0 otherwise a longer data record would be required. 

Figure 6 shows the identified parameter convergence for y * .2 and low noise 
for the periodic coefficients 0^, © 2 and O^. Excellent identification is shown 
with convergence occurring less than 1 rotor revolution. The computed standard 
deviation bands are also shown in fig. 6 which indicate excellent accuracy. 

Identification results at y » .5 . - Fig. 7 shows the A-matrix elements 
and a 22 for y * .5. The periodicity is considerably greater than for y ■ .2 and 
a finite but still very small reverse flow effect is present but not very 
noticeable in fig. 7. The reverse flow is in effect for ift > 180° (1/2 rotor rev) 
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and results In approximately 10X change in magnitude of the a. v> coefficient at 

*. *• 

3/4 of a rotor revolution. This influence Is still too small to he Identified 
accurately over 1 or 2 rotor revolutions in the presence of measurement noise. 

As a result the reverse flow coefficients iL and were ignored in the 
Identification and only 0., 0,, and CL are Identified. 

Figure 8 shows the Identified parameter convergence for the periodic 
coefficients 0^, 0.,, and 0.^. Identification accuracy is excellent and convergence 
to simulation values occurs within l rotor revolution for the low measurement 
noise case. Initial parameter estimates are set to aero with large initial 
covariance. Convergence is rapid and excellent for all three parameter values. 

Ident ifica ti on Results at p « .8, inc ludi ng revers e flow. - The coefficients 
of t he A-matrix a.,^ and a„., are plotted in figure 9 for advance ratio u • .8. At 
this condition very strong periodicity exists including reverse flow effects. 
Figure 9 shows the region of reverse flow and its effect is to flatten the sine 
wave response for 180° v v «. 3h0° (1/2 rev. to l rev. of the rotor). 

The parameter vector used in the identification includes all the periodic 
and reverse flow coefficients of the A-matrix elements a. and a.,., of Kq. (Lh) 
and Kq. (J.7). The parameter vector is 

0 T - [0j CL 0., 0, 0 5 1 (3.12) 

whore the constant terms a,^ and a. v> are assumed known from the hover identifica- 
tion results. The Identified parameter convergence Is shown in figure 10 for the 
sin^> % cos^ terms* figure 11 for the sin^ 1 cosil 1 term and figure 12 for the reverse 

flow coefficients 0, and 0 C . The periodic coefficients 0, * 0* and shown in 

s > I 2 J 

figure 10 and figure 11 show excellent convergence as tuned by the standard 
deviation band + lo. The reverse flow terms 0^ and 0^ shown in figure 12 show 
s’ow convergence as noted hv the + lo bands. The slow convergence is due to the 
low information content of the reverse flow coefficients which is in effect over 
only 1/2 of a rotor revolution. Improved accuracy can be obtained by either 
Increasing the signal to noise ratio or averaging over manv revolutions. 

N V 

Figure 1J shows the parameter convergence for 0^ and 0 ^ for an increase 
in the signal level by a factor of J. Convergence is better than in figure 12 
as shown bv the smaller + lo bands. Parameter convergence is further shown in 
figure 14 for the increased signal level (xJ) for a second rotor revolution 
using the estimated parameters and covariance at the end of the first revolution 
as initial estimates. Convergence is excellent after 2 rotor revolutions as shown 
In figure 14. 
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The results at advance ratio y ■ .8 have accurately identified the periodic 

coefficient model of Eq. (3.1) through Eq. (3.3), including the reverse flow 

terms. Since the identification is from an initial condition response the control 

coefficients of Mq shown in Eq. (3.4) are not identified. The next section 
R 

discusses the results of the stochastic control application which identifies both 
the A-matrix elements and control matrix elements while exercising control. 

Stochastic Control Results 

The stochastic control formulation was briefly outlined in section II which 
showed the objective function for minimization (Eq. (2.3)) and the dynamic model 
state space description in Eq. (2.1) and Eq. (2.2). The detailed discussion of 
the control of systems under uncertainty, the certainty equivalence (CE) property, 
the caution property, and the probing property was also presented. Based upon 
these concepts, two cases are considered. First, the CE control is used at 
U = 0 and the parameter identification convergence and control performance is 
discussed. Second, the periodic A-matrix elements are treated as a random walk 
model and an assessment of the caution and probing significance is addressed. 


Stochastic Control results at y « 0 using the CE controller . The CE control 
algorithm is shown in Eq. (2.57) through Eq. (2.59) where the parameters 0(k) are 
identified on-line with the EKF. Since accurate parameter convergence requires 
1/4 to 1/2 of a rotor revolution in time (.2 seconds), the control solution (Eq. 
(2.57)) is expected to perform less satisfactorily in the beginning and then 
improve as time goes on. The simulation model is shown in Eq. (3.5) where the 
parameter vector for identification Is 


0 


l*. 


‘22 


b 21> 


(3.13) 


and an initial condition 8 ■ 35 deg/sec is used on the system. 

The control objective is to minimize the criterion 
/ N 

J * e ]T (q^ x^(k) + q 2 x^(k) + ^ u^(k)) + q N x^(N) + q N x^(N) 
' k-1 12 


(3.14) 


where, N * 41, q^ * 0, q 2 = 0, * .2, q N * 1.0, and q^ - .01. 

The weighting terms in the performance criterion were determined by trial and 
error such that the feedback was not so large as to prohibit accurate parameter 



J 


identification. The state weightings and q£ were set to zero to permit state 
response of large enough magnitude for accura e parameter identification. The 
optimization was performed over N * 41 time steps which corresponds to 1/2 rotor 
revolution. The control weighting r^ was selected to keep the control magnitude 
reasonably small. 

Figure 15 shows the identified parameter convergence during the CE feedback 
control operation. The parameters are identified on-line and the controlled 
state response is shown in figure 16. The uncontrolled response (IC response) 
and the controlled response are compared in figure 16. The stability of the 
controlled system is increased as shown by the faster decay rate. This increased 
stability is further shown in figure 17, '.ch slows the equivalent closed loop 
system damping parameter 3 22 q L ^ a 22c L = a °2 “ ^21 L 2^ ’ t * me var yi n 8 feedback 

gain increases significantly near the final time causing fu-r eased stability 
as time approachs N = 41. Note that this was dictated by selecting zero weights 


(q^ * 0, q£ ■ 0) on the state and nonzero weights (q N ^ 


1., q„ * .71) at the 
N2 


terminal time. This selection permits better parameter identification to occur 
since the feedback is small initially and the state response is nearly that of the 
uncontrolled system. A time response of the control u^ is shown in figure 18 
which shows increasing control magnitude near the terminal time. 

In summary, the CE control was found to accurately identify the parameters 
82^ and 822 while providing increased stability over the uncontrolled system. 

The parameter b2^ requires further data length to improve upon accuracy. The 
feedback control gains were computed based upon current identified parameter 
estimates. The CE control was found to perform in a satisfactory manner, that 
is it provided increased stability for an uncertain system. The parameters of 
the system were identified while being controlled. Passively adaptive control 
proved successful for the constant coefficient system at hover. 

To examine the contribution of the deterministic, caution and probing 
control, the performance criterion is decomposed as computed using the concepts 
presented in section II and ref. 26. The cost decomposition for the hover case 
when the model is assumed to have constant parameters is 


J “ J DET + J PR0B + J CAUT 


(3.15) 


where , 


J DET " 
J PR0B = 
J CAUT = 


.796 

.012 

.181 


J * .989 
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As shown from the cost decomposition the probing component is negligible and 
caution accounts for less than 20% to the total cost J. Therefore, the control 
solution is dominated by the deterministic aspects of the problem and this 
explains the excellent performance of the CE control (probing and caution aspects 
are not included in the CE control design) , 

Stochastic control results for y > 0 . - The dual control solution was 
investigated for this case by examination of the cost decomposition as was done 
at u * 0. A random walk model was assumed for the parameters to account for the 
periodic variation of the coefficients a^, and b^. The results of the 
cost decomposition for the case y * .8 are shown in Eq. (3.16) 


where, 


J * J DET + J PR0B + J fAUT 


J DET 

J PR0B 

J CAUT 

J 


1.92 

.582 

3.28 

5.79 


(3.16) 


The caution component dominates the performance and is approximately 60% of 
the total cost. The deterministic component is nearly 30% and probing accounts 
for 10% of the total. These results imply that the control solution should include 
the caution property (the probing property can reduce the cost further, but not 
more than 10%) . 

Comparing the cost decomposition of y * 0 with y > 0 indicates that the 
time variation of parameters requires caution and probing properties. The 
larger the variation in the parameters the more dominate would be the caution 
and probing properties for successful identification and control. 
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GROUND RESONANCE PARAMETER IDENTIFICATION 


The results of application of parameter identification to helicopter ground 
resonance is presented in this section. The EKF is used as the parameter identi- 
fication method and results of a multiblade coordinate measurement model is 
assumed as well as a single blade measurement. The single blade measurement 
results in periodic coefficients in the measurement equations. 


Ground Resonance Model Description 


The parameter identification objective is to identify the rotor blade damping 
and thus, using a linearized multiblade coordinate model, the system modal damping 
(eigenvalues) can be determined. This differs from the modal damping determination 
technique of ref. 35 which determines the modal damping directly using the moving- 
block Fast Fourier Transform (FFT) approach. In ref. 7 modal damping is deter- 
mined from random responses. The advantage of the state space approach is that 
the blade damping, can be identified at a lower RPM condition (i.e. before 
coalescence with support frequency). Then, stability (eigenvalues) can be pre- 
dicted at all RPM conditions before testing at these conditions. This capability 
permits considerable safety in testing new rotor designs. The modal damping 
identification techniques of ref. 35 and ref. 7 can only be applied at the RPM 
condition under test (i.e. prediction capability is not possible). 

The second order equations of motion for the ground resonance description 
include multiblade rotor coordinates C^ c and for rotor lag degree of freedom 
and the support x and y degrees of freedom. The model is shown in Eq. (3.17) 


and is representative of a full scale rotor system being tested on a wind tunnel 
support system. 
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These equations result from the original derivation by Coleman [36] and are also 
developed in ref. 37. Ground resonance is of concern when the rotor RPM, is 
such that the inplane regressing frequency coalesces with the support degrees of 
freedom. This can result in an instability and is thus a concern during heli- 
copter wind tunnel tests of new rotor designs. 

Two measurement systems are addressed. The multiblade measurement model 
assumes measurement of 4c’ 4s’ x ’ y (i.e. each degree of freedom). The single 
blade measurement model is 


5 * cos<J> + i^ s sin^i (3.18) 

and also includes the support positions x and y. 

The nominal parameter -ralues for support and rotor degrees of freedom are 


w « 11.6 rad/sec u)_ * .3 

X C, 

W “ 14.6 rad/sec n * .03 

y C 

n * .04 
x 

n = .04 

y 


(3.19) 


The 
model is 


continuous in time state variable description of the ground resonance 
shown in Eq. (3.20; and Eq. (3.21). 


x * Ax + Bu 
y * Cx 

The state, control and nonrotating measurements are defined as 


(3.20) 

(3.21) 


i 

A 



X * 



X 


y 



• (.« * 
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(3.22) 


32 



The rotating measurements and support coordinates are 
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( 3 . 23 ) 


The A and B matrices can be expressed in terms of RPM (ft) and inplane 
damping n g as follows. 


H 
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B - 


0 

0 

0 

0 


2.7146 x10 

U 2 


-5 


0 
0 
0 
0 

3.0570x10 


-5 


or 


1.8097x10 


-5 


ft 


2 


2.0381x10 


-5 


n 
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(3.25) 


The nonrotating measurement matrix is 


C - 


10000000 

01000000 

00100000 

00010000 


(3.26) 


The single blade rotating measurement and support measurement matrix is 


C 


COS 4/ 
0 
0 


sin 

0 

0 


0 0 
1 0 
0 1 


0 0 
0 0 
0 0 


0 0 
0 0 
0 0 


(3.27) 


The state equations defined in Eq. (3.20) through (3.27) are used in discrete form 
in the identification. 


The eigenvalues for four selected conditions are obtained from the A matrix 
of Eq. (3.24). The four conditions yield either a stable or unstable condition 
depending on the value of rotor rotational speed CD and blade inplane damping 
(H J . The four conditions selected are; 


N 


35 


i 


CASE 

ft (RPM) 


COMMENTS 
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1 

150 

.03 

Stable 

2 

150 

.02 

Stable 


3 

160 

.03 

Stable 


4 

160 

.02 

Unstable 



The complete eigenvalues are shown in Table 3. 

For ground resonance detection, It is desirable to identify and monitor 
continuously the stability as RPM is Increased. Cases 1 and 3 represents a stable 
situation and cases 2 and 4 represents a change from a stable to an unstable 
situation. The stable-to-unstable situation could be predicted without testing 
at the unstable condition if blade damping was identified at case 2 (ft - 150 RPM) 
and then eigenvalues obtained from the A-matrix (Eq. 3.24) for ft - 160 RPM. This 
prediction capability is a primary advantage of the state space approach over the 
modal approach of ref. 7 and ref. 35. 

Use of Free Response Data for Identification of Damping for Ground Resonance . - 
The original goal of this research activity was to utilize the dual property of 
stochastic control to enhance identification of parameters. The emphasis of the 
ground resonance study was toward the identification aspects. Simulations were 
performed using the dual control solution and it was concluded that the dual 
control solution required excessive control magnitudes for the system modeled by 
Eq. (3.17) and E H . (3.19). The magnitude of the available control is found to be 
too constrained in magnitude so that the dual control cannot be effective in pro- 
ducing a significant response over small time Intervals. 

Since, the main emphasis of the ground resonance study is on the parameter 
Identification aspects, it was decided that identification from free response data 
or open loop forced response data would be a more effective solution for the ground 
resonance problem. Thus, all results presented are from free response data for the 
ground resonance problem. 

Numerous simulations were performed using the ground resonance math model. 

A sample time of At - .005 was found to yield accurate results and thus all simula- 
tions were performed with this sample time. In addition, all simulation runs were 
performed at a RPM of ft » 150 and blade in-plane damping of n, * .03. This 
represents a lightly damped case and is representative of a typical condition for 
which parameter identification is important. 
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All identification runs were performed from simulated initial condition free 
response data. Five different initial condition sets were investigated shown in 
Table 4. 

Initial condition set 1 was an arbitrary guess, whereas initial condition set 
2 through set 4 was obtained by simulation of the ground resonance model with sine 
input forcing on F . Initial condition set 2 through set 4 was obtained from the 
forced state response at different time points. This procedure closely resembles 
that used to generate free response data using stick stirring for flight and wind 
tunnel testing. Initial condition set 5 was obtained by applying an impulse to the 
F x control input and computing the state variable response at t • 0+. Figure 19 
through figure 24 shows the state response x and y and measurements C for 

IC set 2, 3 and 5. All figures show 150 time steps (.75 secs). 

Random measurement noise was added to all measurements and three different 
noise sequence sets were investigated. The noise levels for the constant and 
periodic model are as follows. 

Standard Deviation Standard Deviation 

of Measurement Noise of Measurement Noise 


(constant 

Coeff. Model) 

Si 

>eriodic 

Coeff, Model) 

4c 

rad 

5.3xl0~ 4 

; 

rad 

5.3xl0~ 4 

4s 

rad 

5.3xl0~ 4 

X 

- * 

3.9xl0 -5 

X 

_ * 

3.9xl0 -5 

y 

_ * 

3.9xl0 -5 

y 

- * 

3 . 9xl0~ 5 





* Nondlmensional (-r R) 

The noise levels were selected such that 30 was approximately 10% of the 

signal level. For example, for an assumed blade lag response of 1 cepcee 

-4 

(.017 rad.), 3 ■ 15.9x10 is approximately 10% of lag response. 

Two different initial starting estimates for the blade lag damping were used. 
A initial estimate f|j- * .04 was selected to be close to the true parameter used 
in the simulation (0j rue ■ .03). An initial damping estimate of hj. * .1 was used 
to repr'-jent an initial guess far removed from the simulation value. 

Parametric studies were performed using the conditions outlined above to 
study parameter identification convergence for blade damping. 
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Parameter Identification Results 

This section shows the results of parameter identification from the ground 
resonance simulation model which was outlined in the last section. All computer 
runs were performed at 150 RPM and blade inplane damping ■ .03. The identifi- 
cation runs were performed using up to a maximum of .75 seconds of data (approxi- 
mately 2 rotor revolutions) at a sample interval of At ■ .005 sec. The lightly 
damped mode has a frequency of 11.07 rad/sec (.707fl) and thus .75 seconds of data 
is approximately 1 1/2 cycles of this mode. 

The results are presented in Figure 25 through Figure 33. Table 5 is presented 
for convenience which summarizes the conditions used for each run presented in the 
figures. Each figure will be discussed separately with reference to Table 5. 

Figure 25 shows identified damping convergence for three differ, .it random 
measurement noise sequences. Initial condition set 1 (the arbitrary guess) ie 
used and the parameter process noise covariance Q 0 ■ 0. The top figure shows the 
results for the three noise sequences for only 40 samples (.2 secs). The bottom 
figure is a continuation for noise sequence s :t 3 and is shown for 80 samples 

A 

(.4 sec). The initial parameter starting value is - .04 und is close to ^he 
true value of ■ -03, Convergence is good for all three noise sets and renuires 
80 samples or more for convergence as shown by the bottom figure. 

Figure 26 shows identified damping for three different initial condition 
sets (Q 0 * 0). Initial condition set 5 has more information content for improved 
identif ibility whereas IC set 2 has less information content thus convergence is 
slower. The rate of convergence is reflected in the covariance of the estimate 
(+ la is shown in the figure). The convergence of IC set 4 (denoted by the circle) 
has falsely converged to .042 because the estimate of the damping has deviated far 
from the true value resulting in violation of linearizations in the EA* IC set 4 
is an example of poor convergence as discussed in section 2.0 because the process 
noise :ovariance Q 0 -■ 0. IC set 2 and IC set 5 show acceptable convergence. 

Figure 27 shows identified parameter convergence comparing the constant and 
the periodic coefficient measurement model. Again the initial damping n *■ .04 

r\ ** 

and Q *0. Convergence is good in both cases with the constant coefficient 
model showing fater convergence. A longer data record is required for the 
periodic model. 
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Figure 28 through Figure 33 use an Initial parameter estimate of * •! 

which is far removed from the simulation value of * .03. The constant 

coefficient model results are shown in Figure 28 (Q® ■ 0) and Figure 29 ( 0 ® f 0) . 

False convergence is shown in Figure 28 with Q® * 0, whereas good convergence is 

shown in Figure 29 with Q® i 0. The + 2a uncertainty bands in the figures reflect 

> 

convergence. 

Figuie 30 and Figure 31 show results for the periodic coefficient model for 

Q n 

Q “0 and Q v ? 0, respectively. The results are similar to the constant coefficient 
case of Figures 28 and 29. Excellent convergence is shown in Figure 31 which uses 
Q® ^ 0. IC set 5 is used with .75 seconds (150 samples of data). 

Figure 32 shows results for the periodic measurement model with initial condi- 
tion set 2 obtained from stick stirring. Convergence is acceptable as determined 
by the + 2a band shown in the figure with Q® 4 0. Convergence to the true parameter 
value requires a longer data record. 

Figure 33 shows results for the periodic measurement model with IC set 3 
also obtained from stick stirring. IC set 3 was selected after allowing a further 
increase in the time over IC set 2 before setting the control input to zero. The 
initial conditions thus obtained yield better information for the resulting free 
response data. Therefore, identified parameter convergence is faster for IC set 
3 over IC set 2. This is reflected in the + 2a band shown in the figures. Con- 
vergence is excellent as shown in Figure 33. 

Figure 33 demonstrates that parameter identification of inplane blade damping 
is possibl using less th?n 1 second of data (approximately 2 rotor revolutions) 
using a single blade measurement of blade angle in the rotating system (periodic 
measurement model) and x and y support position measurements. 

A summary of main conclusions for the ground resonance study is as follows; 

. The use of free response data for identification of the in-plane damping 
parameter has been shown to yield excellent convergence using stick 
stirring to obtain initial conditions. 

. Excellent parameter identification of damping was shown for both the 
constant and periodic measurement model. 

. The extended Kalman filter (for the no process noise case) has shown 
excellent identified parameter convergence (providing the process noise 
covariance on the parameter, Q®, is set to a nonzero val_:e) . 

. The state variable multiblade model A-matrix can be used to predict the 
stability for RPM conditions before they are tested. This feature makes 
this approach to identification very competative with the modal approach 
of ref. 35 and ref. 7. 
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CONCLUSIONS 


A general procedure for parameter identification and stochastic control of 
helicopter rotor dynamic systems has been presented. The formulation is based 
upon a state variable representation described by constant and periodic coeffi- 
cients. The extended Kalman filter was shown to yield unbiased parameter estimates 
for all case studies for both periodic and constant coefficient models. This was 
found to be true for systems without process noise. 

The pure identification results of the rotor flap model were successful over 
an advance ratio of y ■ 0 to (i ■ .8. This included constant and periodic coeffi- 
cients and the effects of reverse flow. Typically less than 1 rotor revolution 
was required where there was sufficient signal to noise ratio. Stochastic control 
of the rotor flap model was successful in that stability was improved while para- 
meters were accurately identified. The CE control was found to be successful at 
■ 0. For systems with time varing coefficients, the caution property was 

found to dominate, with the probing property contributing 10% to the total cost 
function. 

Pure identification was successful for the ground resonance problem for 
either the multiblade coordinate rotor inplane measurement model or for a single 
rotating blade measurement. Typically 2 rotor revolutions are required in 10% 
measurement noise where free-response data is used (generated from a stick 
stirring control input). The state variable identification approach has the 
advantage over modal methods (e.g. the moving-block FFT method) in that stability 
prediction over the complete RPM range is possible from test data obtained at one 
test RPM condition. 

The EKF identification method is an efficient identification approach for 
systems without process noise. The "corrected" EKF (Ljung 1979) should be 
investigated for systems with both measurement and process noise. The corrected 
EKF method could oe assessed via simulation and test data applications. 

The rotor flap model Investigated iu this study could be extended to include 
blade inflow models and higher degrees of blade flexibility. Identification of 
periodic coefficient models should include a hypothesis testing procedure (e.g. 
stepwise multiple regression) for automatic determination of the required periodic 
terms in the model. This is important since, ir. general, it is not known which 
periodic terms would be required, particularly in the reverse flow region. The 
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identification and stochastic control methods u J in this research should be 
investigated as to more efficient computational algorithms, which is of particula 
importance for real time or near real time applications. 
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APPENDIX A 

AUTOMATIC CALCULATION OF THE JACOBIAN 

The Extended Kalman Filter is quite similar to the standard Kalman Filter 
with the exception of the need to evaluate the Jacobian (2.19) at every step. 

In practice the unknown parameters appear in several elements of the contin- 
uous time system matrices A(0) and B(0) from (2.4). Furthermore, a particular 
parameter usually appears in se\eral elements of the system matrices. 

The following procedure, applicable for a wide class of problems, can be 
used to automate the calculation of the Jacobian. 

It is assumed that the system matrices A and B from (2.4) have elements a^ 
and by that are either known or given by the product of a known constant with an 
unknown parameter. To indicate which parameters enter where in A and B define 
the matrices 



(A. 1) 

B p - 

(A. 2) 

of dimensions n x xn x and n x xn u > respectively, where 


( 0 if a . , = a , . 

(A. 3) 

( m if a. J 0 m 



where 1 < m < n^ and 


* ■ u u 1 


(A.4) 


is a known matrix and similarly for . . 
Thus 






a. .0 
ij m 


if a*r . = 0 
ij 

it m « a^j # 0 


(A. 5) 


Propagation of the state (2,15) requires the evaluation of 
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f [z(k|k) , u(k) ] 


OFpqqd is 

F(0(k|k))ic(k|k) + 9^^)u(k)l 


0(k|k) 


(A. 6) 


which, in turn, requires the evaluation of F(0(k|k)). The latter can be then 
easily obtained based on (A. 5) as 


N 


F(0(k|k)) = I + A(0(k|k))A 


(A. 7) 


with 


a (0(k|k)) 




i ij 0 m (k|k ) 


if a 


if 


ij 


H 


m * a; . / 0 


(A. 8) 


The covariance propagation equation (2.18) requires the evaluation of the 

Jacobian (2.19) which is an n xn matrix. Based on (A.1)-(A.4) the expressions 

z z 

of the elements of the Jacobian matrix 


*(k) * [*. 4 (k)] 
ij 


(A. 9) 


are 


$^(k) * Fjj (0(k k)) i=l,...,n x ; j-l,...,n x 


$ ij (k) " 1 


i=j=n x + m ; m«l,...,n fl 


0 


(A. 10) 
(A. 11) 


I T 

- - - (k) * \.rt.. . -a 


i,n + m 
x 


+ z 


'S.eL (i,m) 


l jn X f (k|k) nr -y » 1A 

JleL (i,m) 


u j(h| 


V k) 


i ® !»•••» » m® 1,..., n Q 


* 0 otherwise 


0 


(A. 12) 
(A. 13) 


where 


L a (i,m) 


u • - p 


a 'u ■ m) 


(A. 14) 
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is the set of state components that multiply 0 in the i-th row of (2.12) and 

in 


L b (i,m) ■ (J, : = m} 


(A. 15) 


is the set of control components that multiply 0^ in the i-th row of (2.12). 

In a similar manner one can compute the Hessian of f^, the 8,-th component 
of the vector valued function f from (2.12), denoted as 


_ 


9 

.£ 

^ Z 9 zz 


dZj 

L. 


. * 

The expression of the elements of this symmetric matrix 

.P 


(A. 16) 


a^ A if a*^ » m for i = l,...,n ; m = l,...,n p 


0 


\jj m \b 

n +m,i i.n +m 
x x 


(A. 17) 


0 Otherwise 
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Table 1. - Identification Methods Which Yield Asymptotically Unbiased Estimates 
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Table 2. - Flap Equation Periodic Coefficient Va’ues with Advance Ratio 
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Reverse Flow 


Table 3. - Eigenvalues for Ground Resonance Model (Frequency Nondimensloned by U) 
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IC set 2 through 4 obtainea - 

Tr cet 5 obtained from pulse input. 



Table 5, — Parameter Identification Run Summary for Ground Resonance Study 
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2 frou stick stirring, proper 
ence but requires more data 
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Figure 1. - Stochastic Controllers Based upon the Certainty Equivelance 
Property and the Dual Property of Stochastic Control Theory. 


IC RESPONSE 




Figure 2. - Rotor Flap M del State Time History Response Obtained From 
Free Response Portion of Pulse Input. (6(0) = 35 deg/sec). 
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Figure 3. - Rotor Flap Model Identified Parameter Convergence, p * 0 
(No Periodicity, High Noise). 
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Figure 4. - Rotor Flap Model Identified Parameter Convergence, p * 0 
(No Periodicity, Low Noise). 
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TIME, SECS 

Figure 5. - Rotor Flap Model Periodic Coefficient Variation of A-Matrix 
(U - .2, Small Periodicity, Negligible Reverse Flow). 
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Figure o. - Rotor Flap Model Identified Parameter Convergence, p ■ .2 
(Reverse Flow Cooff. Not Identified). 



1 

A 


1 

! 

| 


ORIGINAL PAGE IS 




0 .1 .2 .3 .4 


TIME, SECS 

Figure 7. - Rotor Flap Model Periodic Coefficient Variation of A-Matrix 
(y * .5, Moderate Periodicity, Small Reverse Flow). 
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Figure 8. - Rotor Flap Model Identified Parameter Convergence, y » .5 
(Reverse Flow Coeff. Not Identified). 
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Figure 9. - Rotor Flap Model Periodic Coefficient Variation of A-Matrix 
(y = .8, Strong Periodicity, Reverse Flow). 



Figure 10. - Rotor Flap Model Identified Parameter Convergence, y = .8 
(cos sin ip Coefficients, 1st Revolution). 
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Figure 11. - Rotor Flap Model Identified Parameter Convergence, y = .8 
(cos ^ sin ip Coefficient) 
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Figure 12. - Ritor Flap Model Identified Parameter Convergence, y * .8 
(Reverse Flow Coefficients) . 
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Figure 13. - Kotor Flap Model Identified Parameter Convergence, y * .8 
(Reverse Flow Coeff., Signal x 3). 
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Figure 14. - Rotor Flap Model Identified Parameter Convergence, y - .8 
(Reverse Flow Coeff., Signal x 3). 
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Figure 15. - Rotor Flap Model Identified Parameter Convergence for the 
CE Stochastic Controller, y - 0. 
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Figure 16. - Rotor Flap Model State Time History Response Comparing 
Open Loop Free Response with CE Controlled P.esponse. 
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Figure 17. 
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Figure 18. 
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Rotor Flap Model Equivalent Closed Loop System Damping, CE 
Control, u ■ 0 (Time Varying, A rT ■ A-Blc) . 



Rotor Flap Model Feedback Control Time History Response, 
CE Control, p - 0. 
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Figure 19. - Ground Resonance Rotor Model State Response, IC Set 2, ft • 150 RPM. 



Figure 20. - Ground Resonance Model Measurement Response, IC set 2, ft ■ 150 RPM 
(Rotating Rotor Measurement) . 


62 




ORIGINAL P/ ~ , c ft , 

0F P °°« gfcfi 


150 


Figure 21. - Ground Resonance Rotor Model State Response, IC Set 3, ft ■ 150 RPM. 



Figure 22. - Ground Resonance Model Measurement Response, IC Set 3, ft » 150 RPM 
(Rotating Rotor Measurement). 
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Figure 25. - Ground Resonance Model Identified Parameter Convergence Compari" 
Three Random Noise Sequences (Q 0 » 0, IC set 1, ft * 150 RPM) 
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Figure 26. - Ground Resonance Model Identified Parameter Convergence Comparing 
Different Initial Condition Sets (Noi e Set 3, Q 0 - 0, fl - 150 RPM) 
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Figure 27. - Ground Resonance Model Identified Parameter Convergence Comparing 

Periodic and Constant Coefficient Model (It Set 5, 0, Q • 15.' RPM) , 
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Figure 28. •* Ground Resonance Model Identified Damping Parameter Convergence 

Using Multiblade Lag Measurement (Constant Coefficient Measurement 
Model, Prr is Noise Covariance of Parameter Q® - 0) . 
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Figure 29. - Ground Resonance Model Identified Damping Parameter Convergence Using 
Multiblade Lag Measurement (Constant Coefficient Measurement Model") 
Process Noise Covariance of Parameter Q® * 4x10“^) . 
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Figure 30. - Ground Resonance Model Identified Damping Parameter Convergence Using 

Blade Lag Measurement (Periodic Measurement Model) Showing False 
° P 

Convergence when Q u = 0. 
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Figure 31. - Ground Resonance Model Identified Damping Parameter Convergence Using 
Blade Lag Measurement (Periodic Measurement Model) (Process Noise 
Covariance of Parameter Q® * 4x10“^) . 



Figure 32. - Ground Resonance Model Identified Damping Parameter Convergence 

(Periodic Model, Q® - 4x] ''“ 6 , IC Set 2 from Stick Stirring, 0, - 150 RPM) . 
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Figure 33. - Ground Resonance Model Identified Damping Parameter Convergence 

(Periodi~ Model. Q® = 4x10”^, IC Set 3 from Stick Stirring, ft = 150 RPM) 



